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We propose a change of style for numerical estimations of physical quantities from measurements 
to inferences. We estimate the most probable quantities for all the parameter region simultaneously 
by using the raw data cooperatively. Estimations with higher precisions are made possible. We 
can obtain a physical quantity as a continuous function, which is differentiated to obtain another 
quantity. We applied the method to the Heisenberg spin-glass model in three dimensions. A dynamic 
correlation-length scaling analysis suggests that the spin-glass and the chiral-glass transitions occur 
at the same temperature with a common exponent u. The value is consistent with the experimental 
results. We found that a size-crossover effect explains a spin-chirality separation problem. 


Introduction- Estimations of physical quantities in nu¬ 
merical simulations are based on equilibrium statistical 
physics [lj. We virtualize a model system in a computer 
and perform independent measurements on the system 
using a definition of a physical quantity. When an evalu¬ 
ation process is complex, both systematic and statistical 
errors are accumulated in the obtained data. We some¬ 
times encounter numerical instabilities, which may affect 
a final physical conclusion. In what follows, we explain 
the situation of interest using a correlation-length esti¬ 
mation. 

An estimation formula for a correlation length, 
£, is given by the second-moment method: £ = 

Vxo/Xk - 1/A; i- Here, xo denotes the susceptibility 
and Xfc its Fourier transform with k as the lowest wave 
number of the system. This expression itself is problem¬ 
atic. Both numerator and denominator of this expression 
approach zero as the system size increases [k —>• 0), where 
this formula becomes exact. We encounter the numeri¬ 
cal instability caused by the expression 0/0. In order to 
avoid this problem, Bellettiet al. Q proposed the reduc¬ 
tion of this instability by estimating £ through the inte¬ 
grals Ik = Ldrr k f(r) and £ is obtained as h/h- Suwa 
and Todo |4j proposed a generalized moment method for 
gap (A ~ l/£) estimation in quantum systems. System¬ 
atic errors and ambiguity caused by using small-Z, data 
are eliminated. 

Recently, big-data handling has become possible due to 
rapid increase in computational power. Data science is 
now one of the most promising fields in science and tech¬ 
nology. As regards its application to physics, the topic 
of Bayesian inference has attracted considerable inter¬ 
est @,[g|. In this context, Harada @] introduced Bayesian 
inference into a parameter estimation of the finite-size 
scaling analysis. 

In this paper, we extend its application to estimations 
of physical quantities. For example, we can obtain an 
analytic expression for an energy out of the discrete raw 
data as the most-probable model function. Then, we 
obtain the specific heat by analytically differentiating it. 
A critical temperture is estimated automatically within 
this procedure. Since directly-observed (raw) data are 


cooperatively utilized in this inference procedure, we can 
reduce numerical errors and avoid numerical instabilities. 

We also discuss in this paper a size crossover effect in 
random systems. Hukushima and Campbell [sj reported 
that there exists a crossover size, L ~ 24, where the finite- 
size effect of the correlation-length ratio, £/L, changes 
its trend from increasing to decreasing in the Ising spin- 
glass model. Similar non-monotonic size dependences 
have been observed in the ±J Heisenberg spin-glass 
model. The chiral-glass susceptibility of sizes smaller 
than L = 39 increases with the system size but that 
of larger sizes decreases @. Size-crossover effects were 
also observed in a random quantum spin chain 0,0. 


Short-range spin correlations exhibit an exponential de¬ 
cay, which suggests that the energy gap is finite; in con¬ 
trast, the long-range ones exhibit an algebraic decay in¬ 
dicating that the energy gap is zero. The size-crossover 
effect may influence the final physical conclusion. We ex¬ 
plain contradictory arguments on a spin-glass transition 
by this effect. 

Method- We explain the method in a two-dimensional 
Ising ferromagnetic model. We performed equilibrium 
simulations and obtained data for energy, £), and the 
magnetization, Mi, at each temperature, T), where i is 
the data index. The linear system size is 999, and it is 
set to 1999 in the vicinity of the transition temperature. 
These data are depicted in Fig. [T] by circle symbols. We 
fit them by the Gaussian kernel regression [7, 12] using 
three variables, Xi, xji, and e, defined as 


Xi = In \Ti — T c \ 
Vi = ln(-Ej) 
e, = {6 Ei/Ei) 2 


(Vi = In |Mj| for M) 

(ej = (SMi/Mi) 2 for M). 


Here, SEi and SMi denote errors for E and M, and 
T c denotes the critical temperature that is to be esti¬ 
mated in the following analysis. We defined a Gaussian 


= el 


0 exp 


(Xj-Xj 


28 ? 


+ 


kernel function as K (xi , Xj 
where 0q,9i, and 02 are hyper parameters. A general¬ 
ized covariance matrix is Ey = eiSij + K(xi,Xj). Then, 
the following log-likelihood function is to be maximized: 


In £ = — i In | 27 t£ | — 




This function is de- 
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FIG. 1: (Color online) Temperature dependences of 

energy(E), magnetization (M), the specific heat(C), and the 
effective critical exponent /5 e g in the two-dimensional Ising 
ferromagnetic model. Error bars are smaller than the sym¬ 
bols and line widths. 


fined independently for both T) > T c and T) < T c , and 
we take a summation of them. The hyper parameters, 
{ 0 o, 0 i, 02 }, are a ls° defined independently for two re¬ 
gions. We searched for seven parameters, one T c and 
two sets of {9q, 6 \, 62 }, that maximizes the log-likelihood 
function by using the downhill simplex method [l3| ■ We 
tried this search for four hundred times by changing the 
initial values of the parameters. We estimated averages 
and error bars for parameters over them. The critical 
temperature is obtained as a parameter that separates 
the data into two regions, where the data are fitted most 
smoothly. It was T c = 2.26914(4) for the E inference 
and T c = 2.26919(4) for the M inference. They are very 
close to the exact value T c = 2.269185314 • • •. Using the 
obtained parameter set, we write a model expression for 
E as 

E(T) = - exp [K{x,Xi)H~^yj], x = In |T - T c |, (1) 

where the summation over i and j are taken. We dif¬ 
ferentiate this function analytically, and we obtain the 
specific heat, C(T) = as a continuous function. The 
inference results for E(T) and C(T) are depicted by lines 
in Fig. Q] We confirmed that the C(T) function is con¬ 
sistent with the exact results. We obtained a function 
for M in the same manner. Since M ~ (T c — T)@, the 
effective /3 is given by /3 e ff = with y = In \M(T)\. A 
critical exponent ^ = 1/8 is a value at T = T c . A critical 
region, where /3 e ff approximately equals to 1/8, is very 
narrow. 

The nonequilibrium relaxation method 0-E3 was 
proposed to treat large systems in a simple and easy 
manner. This approach has been applied successfully in 
random systems [9l4ll[ 17 -l20l|. The dynamic correlation- 


length scaling method [21] was proposed as a variation 
of this method. We use this method together with the 
inference method to clarify the spin-chirality problem in 
the Heisenberg spin-glass model in three dimensions. 

Model- A spin glass is a disordered magnet charac¬ 
terized by frustration and randomness 22|, |23|. One of 


the most important and unsolved problems in spin-glass 
studies is the coupling or separation of the spin-glass (sg) 
degrees of freedom and chiral-glass (eg) degrees of free¬ 
dom [lil I24l - l4ll ] . Kawamura 25, 2fj] introduced the chi¬ 
rality scenario, wherein the eg order exists without the 
sg order. There is another scenario, in which the sg and 
eg transitions occur simultaneously. In 2009, two stud¬ 
ies fii - Eot on this topic drew two opposite conclusions 
even though the authors in each case performed similar 
amounts of simulations, but treated the finite-size effects 
differently. The present situation suggests that we need 
considerably larger system sizes to address this problem. 

Our model Hamiltonian is: H = — y~/ ? - rj) JijSiSj. The 
summation runs over all the nearest-neighbor spin pairs. 
The interactions Jtj take on two values, ±J, with the 
same probability. The temperature T is scaled by J. 
The model is defined on a simple cubic lattice of the 
form N = LxLx(L+l) with L = 159. The skewed 
periodic boundary conditions were applied. We calcu¬ 
lated the sg/cg susceptibility, Xsg/Xcg , sg/cg correlation 
functions, / sg //cg, and sg/cg correlation length, £ S g/£ C g- 
One Monte Carlo (MC) step consists of one heat-bath up¬ 
date, 1/20 Metropolis updates (once every 20 steps), and 
124 over-relaxation updates. All the random bond con¬ 
figurations are different at each temperature. A typical 
sample number at one temperature is 20. More samples 
are treated near and above the transition temperature. 
In the study, we ran simulations at 42 sets of tempera¬ 
tures, and the total sample numbers were 1168. We eval¬ 
uated the order parameters using 435 overlaps among 
30 real replicas. At some lower temperatures, we eval¬ 
uated them using 1128 overlaps among 48 real replicas 
and checked for consistency regarding the replica num¬ 
ber. In the nonequilibrium relaxation study on the spin 
glasses, the thermal average is replaced by the replica av¬ 
erages 21]. The replica number needs to be larger than 
the value in the equilibrium simulations. Numerical error 
bars were estimated in regard to the sample average. 

Results- We prepared the relaxation data of correlation 
functions, f(r,t ) [r denotes distance and t the measur¬ 
ing time step], obtained in the conventional measurement 
scheme. Figure [2] shows the sg and eg correlation func¬ 
tions for typical time steps in the range from t = 316 to 
79433. The temperature, T = 0.200, is close to the tran¬ 
sition temperature. We also plot the small -L (L = 39) 
data at T = 0.210 as shown by circles. The inverse of the 
slope of the curve in this figure corresponds to the cor¬ 
relation length. Here, we found the crossover distance, 
r c , dividing the short-range correlation region and the 
long-range correlation region. Short-range correlations 
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FIG. 2: (Color online) Correlation function data for selected 
time steps for L = 39 (circles) and L = 159 (lines). Arrows 
depict crossover distances between short-range correlations 
and long-range correlations (r c = 9 for sg and r c = 3 for eg). 



FIG. 3: (Color online) Scaled correlation function data for 
r < L/3. 

do not depend on f, T, and L. Meaningful information 
is not included in this region. The growth of the cor¬ 
relation length is only reflected in the long-range corre¬ 
lations. The sg crossover distance (r c ~ 9) is roughly 
three times greater than the eg one (r c — 3). The ef¬ 
fects of the periodic boundary conditions appear as the 
distance approaches L/2. We use only the data in the 
distance range of 2 r c < r < L /3 to exclude influences of 
short-range correlations and the boundary effects. 

The correlation lengths are related to the correlation 
functions via the following scaling ansatz: 

=Hr/m), ( 2 ) 

where T denotes the scaling function and ry e ff is the ef¬ 
fective scaling exponent. In a Gaussian kernel regression 
procedure, we set = r/£(t), yi = f(r, t)/[£(t)] _1_I ' eff , 
and ei = 5yf, with i denoting an index number for all 
the combinations of (r,f). We estimate £(f) and r] e g as 
parameters so that all the f(r,t) data fall onto a single 
scaling function T. Dozens of £ (t) data sets are obtained 
simultaneously from thousands of f(r,t) data sets. Con¬ 
sistency among many data sets yields accurate estimates 


FIG. 4: (Color online) The correlation length data esti¬ 

mated by the inference and that by the measurements (the 
second-moment method). Error bars for the inference data 
are smaller than line widths. Inset: Temperature depen¬ 
dences of the effective exponent r/ e ff. The critical exponent 
77 obtained by the dynamic correlation-length scaling analysis 
is also plotted. 

of the correlation length. 

Figure [3] shows the result of scaling. We rescaled 
£(t) so that the slope of this plot becomes unity as 
T{r/£,{t)) ~ exp[—r/£(f)]. This rescaling defines the unit 
of the length scale. Figure [4] shows the obtained £(t). 
We compare our inference results with those obtained in 
the measurement sheme (the second-moment method). 
The sg data obtained with both methods show a close 
consistency. On the other hand, numerical instabilities 
are observed in the eg estimations by the measurement 
method. In contrast, the inference method solves this in¬ 
stability problem. The effective exponent, rj e g, depends 
on the temperature reflecting a correction to scaling. We 
plot the r/ e ff values in the inset of this figure. It coincides 
with the critical exponent at the transition temperature, 
which will be obtained by the scaling analysis. 

W e ap ply the dynamic correlation-length scaling anal¬ 
ysis |21| using the obtained (£(t),x(t)) data sets. Figure 
[5] shows the scaling plot of the sg and eg transitions. 
We applied the /3-scaling method proposed by Campbell 
et al. [42|. We estimated the scaling parameters by the 
Bayesian inference introduced by Harada j7]. There are 
1187 data sets in this figure, and we chose 800 data sets 
randomly and estimated the scaling parameters for 100 
times. We determined the average and the error bar over 
them. We also plot in the inset the scaling result us¬ 
ing £ obtained by the measurement method. While it 
is impossible to perform scaling analysis on the eg data 
in the measurement method, our inference method made 
it possible. Estimated transition temperatures and crit¬ 
ical exponents are summarized in Table [I] The critical 
temperatures are consistent with previous results. Our 
values of u sg and r] s g are also consistent with those of 
the canonical sg materials [43l - [45| . This evidence sug¬ 
gests that the Heisenberg spin-glass model explains the 
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FIG. 5: (Color online) A scaling plot of (a) the sg transi¬ 
tion and (b) the eg transition. Data of 42 sets of temperature 
ranging from 0.170 to 0.220 are depicted with different color 
lines. Insets depict the same scaling plot using the correla¬ 
tion length data obtained by the measurements (the second- 
moment method). 



T sg 

T cg 

Vsg 

V c g 

Vsg 

Vcg 

this work 

0.203(1) 

0.201(1) 

1.49(3) 

1.53(3) 

0.28(1) 

0.66(1) 

Ref- [18] 

0.21(1) 

0.22(1) 

1.1(2) 


0.27 


Ref- [35] 

0 

0.19(1) 


1.3(2) 


0.8(2) 

Ref-[9] 

0.203(1) 

0.200(1) 

1.79(2) 

1.57(3) 

0.19(1) 

0.83(2) 

Ref- [43] 



1.40(16) 


0.46(10) 


Ref-[44] 



1.30(15) 


0.4(1) 



TABLE I: Comparisons of our results with previous estimates. 
Refs- (is] and -]44| are experimental results of AgMn. 


experiments. 

Discussion and Summary- The evaluations of physi¬ 
cal quantities in numerical studies are generalized to an 
inference scheme. This is a change of style in numer¬ 
ical investigations on statistical physics. We obtain the 
most-probable expression for a physical quantity from the 
discrete raw data. Then, we differentiate or integrate it 
analytically or numerically to obtain various quantities. 
We can improve accuracies of physical quantities because 
they are the product of consistency among many raw 
data sets. This method has potential applications not 
only to numerical studies on theoretical models but also 
to analyses on experimental data. 


In our study on a Heisenberg sg model, we observed a 
simultaneous sg- and cg-transition with a common value 
of exponent v. Here, one may ask why the sg- and cg- 
transitions have been observed sometimes differently and 
sometimes simultaneously. In what follows, we clarify 
this point. There are two important length scales when 
we discuss the finite-size effect. One is the correlation 
length and the other one is the crossover length. In the 
ferromagnetic Heisenberg model, the crossover length is 
only 2-3 lattice spacings. Thus, finite-size scaling analysis 
is possible using data with the minimum size L = 6 [461 ] 
or 8 [47]. As shown in Fig. [3 the sg crossover length is 
9-10 lattice spacings in the Heisenberg spin-glass model. 
This value is comparable with the correlation length in 
the present simulation. The necessary length scale should 
be doubled or tripled under the periodic boundary con¬ 
ditions. This corresponds to a minimum lattice size 
L = 20 — 30. However, these sizes have been mostly the 
maximum sizes in the equilibrium simulations. On the 
other hand, the eg crossover length (r c = 2 — 3) is almost 
same as that in the ferromagnetic model. The necessary 
size may be L = 6 — 8, which has been considered in 
the equilibrium simulations. This crossover-length issue 
is the reason why the sg transition was not detected in 
early simulations, while the eg transition was easily de¬ 
tected. 
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